The present invention relates to biomagnetic field measuring method and apparatus for measuring a biomagnetic field generated by a nerve action of brain as well as a myocardial action of heart of a living body by using a plurality of fluxmeters each consisting of a highly sensitive superconducting quantum interference device (SQUID).
In addition to a magnetic field generated by a current dipole, a magnetic field due to a volume current flowing in the living body is enumerated as a biomagnetic field. Measurement of a normal component (Bz:z component in the Cartesian coordinate system or Br: radius component in the polar coordinate system) is considered to be hardly affected by the volume current. In conventional techniques, the plane of a detection coil connected to a SQUID is disposed in parallel to the body surface to measure Bz or Br which is a normal component vertical to the body surface. Results of the biomagnetic field measurement are displayed in the form of a temporal change waveform of the measured field component or an isomagnetic field map (contour map) for connecting points at which magnitudes of the magnetic field component measured at desired time points are equal to each other. Various analysis methods have been proposed which analyze a magnetic field source participating in generation of the biomagnetic field from the obtained isomagnetic field curve and in typical one of them, analysis is carried out by replacing the magnetic field source with a current dipole.
An isomagnetic field map of a normal component (Bz or Br) of the magnetic field generated by a current dipole is of a pattern having a source pole of magnetic field and a sink pole of magnetic field at positions which are separate from each other from the center where a magnetic field source (current dipole) is positioned. The magnitude, position and direction of the magnetic field source (current dipole) are analyzed in accordance with magnitudes of magnetic field at the two poles and a distance therebetween.
In a first prior art (H. Hosaka and D. Cohen: J. Electrocardiology, 9 (4), pp. 426-432 (1976)), a method is employed for displaying current sources distributed in the myocardium by using an isomagnetic field map of measured normal component Bz with the aim of promoting visibleness of direction and intensity of currents in the myocardium and according to this method, an arrow map is contrived for expressing a current vector J (x, y) defined by equation (1) on measuring points by using an arrow. In the following description, Gothic characters are used to indicate vectors.
xe2x80x83J(x, y) =(∂Bz(x, y)/∂y)exxe2x88x92(∂Bz(x, y)/∂x)eyxe2x80x83xe2x80x83(1)
In equation (1), ex designates a unit vector in x direction and ey designates a unit vector in y direction. This prior art, however, encounters a problem that when a plurality of current sources exist, it is difficult to discriminate the individual current sources from each other on the basis of the isomagnetic field map of normal component Bz.
In a second prior art (K. Tukada et al: Review of the Scientific Instruments, 66(10), pp. 5085-5091 (1995)), for the sake of visualizing a plurality of distributed current sources, the normal component (Bz or Br) is not detected but tangential components Bx and By are measured by using a detection coil whose plane is disposed vertically to the body surface. Each of the measured tangential components Bx and By is displayed in the form of an isomagnetic field map. The tangential components Bx and By measured according to the second prior art are considered to be affected by the volume current but in an isomagnetic field map of two-dimensional vector magnitude Bxy obtained by synthesizing Bx and By measured at time point t pursuant equation (2), a peak can always be obtained directly above a current dipole and therefore, even when a plurality of current dipoles exist, individual current dipoles can be separated for visualization.
xe2x80x83|Bxy(x, y, t)|={square root over ( )}{(Bx(x, y, t))2+(By(x, y, t))2}xe2x80x83xe2x80x83(2)
In a third prior art (Y. Yoshida et al: Tenth International Conference on Biomagnetism, Santana Fe, N. Mex., Feb. 17 (1996)), a normal component and two tangential components of a biomagnetic field are detected by using a vector magnetic field sensor consisting of three detection coils having coil planes which are orthogonal to each other, detection results of the magnetic field components are converted in terms of the Cartesian coordinate system to determine Cartesian coordinate system components Bx, By and Bz, and an isomagnetic field map of the normal component Bz and an isomagnetic field map of two-dimensional vector magnitude Bxy are displayed, respectively.
In a fourth prior art (K. Tsukada et al: Tenth International Conference on Biomagnetism, Santana Fe, N. Mex., Feb. 17 (1996)), two tangential components Bx and By of a biomagnetic field are detected and an isomagnetic field map based on |Bxy|=|Bx+By| is compared with an isomagnetic field map based on a normal component Bz.
Available as diagrams for indicating measurement results of electrical physiological phenomena in a living body are a magnetoencephalogram (MEG) obtained through measurement using a magnetoencephalogram and an electrocardiogram (ECG) obtained through measurement using an electrocardiograph. In measurements of the electrocardiogram, a body surface potential map for mapping an electrocardiographic figure by using a plurality of electrodes is of a well-known technique. The MEG or the body surface potential map is depicted in the form of an isopotential map for connecting isopotential points.
In a fifth prior art (T. J. Montague et al: Circulation 63, No. 5, pp.1166-1172 (1981)), an isointegral map obtained by integrating a temporal change waveform of an output of each one of a plurality of electrodes over a desired time interval is depicted as a body surface potential map.
In the following description, xe2x80x9cbiomagnetic fieldxe2x80x9d means xe2x80x9cmagnetic field generated from a living bodyxe2x80x9d, xe2x80x9ccardiac magnetic field measurementxe2x80x9d means xe2x80x9cmeasurement of a magnetic field generated from the heartxe2x80x9d, and xe2x80x9ccardiac magnetic waveformxe2x80x9d means xe2x80x9cwaveform indicated by a magnetocardiogram (MCG) obtained through cardiac magnetic field measurementxe2x80x9d. Further, xe2x80x9cencephalic magnetic field measurementxe2x80x9d means xe2x80x9cmeasurement of a magnetic field generated from the brainxe2x80x9d and xe2x80x9cencephalic magnetic waveformxe2x80x9d means xe2x80x9cwaveform indicated by a magnetoencephalogram (MEG) obtained through encephalic magnetic field measurementxe2x80x9d.
Each of the conventional isomagnetic field maps of the respective components has inherent features. In the presence of a single current dipole, the position, magnitude and direction of a current source can be analyzed with ease by using the isomagnetic field map of normal component Bz. On the other hand, the isomagnetic field map of two-dimensional vector magnitude Bxy obtained from measurement results of tangential components Bx and By features that even in the presence of a plurality of current dipoles, individual current dipoles can easily be discriminated from each other. But, for detection of a magnetic field, coils are required to be provided in x and y directions and the number of coils is doubled as compared to detection of only the normal component Bz. In vector measurement for measuring all the components Bx, By and Bz, the number of required coils is tripled as compared to detection of only the normal component Bz. Accordingly, the magnetic field sensor consisting of a detection coil and a SQUID is increased in number and in addition, the signal processing circuit and the like are also increased in number, raising a problem that the biomagnetic field measuring system becomes an expensive one. Further, the first prior art is disadvantageous in that arrows are merely indicated on measuring points and detailed distribution states of current sources are hardly discriminated.
From the isomagnetic field map indicated in terms of a biomagnetic field component, the position, magnitude and direction, at a desired time point, of a current source in a living body can be analyzed and detailed information about changes in position, magnitude and direction of the current source can be known. Conventionally, dynamic changes in various kinds of information pieces are captured by using many figures displayed on or delivered to the apparatus so as to diagnose a disease. In the prior arts, however, many diagrams or maps indicating various kinds of information pieces are needed for diagnosis and abnormality of changes in various kinds of information pieces is known empirically. As will be seen from the above, in the prior arts, the processing of displaying, on a single map, systematic information as to what magnitude of current flows through which portion of a living body and as to which region an abnormal bio-current passes through is not executed. In the case of the body surface potential map, an isointegral technique was reported. This isointegral map was drawn by connecting between the same integral values over a desired time interval (for example, an interval during which waves of Q, R and S are generated and an interval during which S to T waves are generated). The advantage of this isointegral map is that information of the heart can be obtained from only a single electrocardiographic figure. But, in the isopotential map when the current source in the heart is assumed to be a single current dipole, a figure results disadvantageously in which an positive peak and a negative peak do not exist immediately above the current dipole but exist at a position which is separate from a point immediately above the current dipole. Further, when the position of the current dipole remains unchanged but the direction of the current dipole changes, the anode and cathode peak positions change, raising a problem that when potential is integrated, correspondence between the current source and the peak of an integral value is impaired. Like the case of the electrocardiogram, mere integration of a component of a biomagnetic field obtained through biomagnetic field measurement faces a problem that the peak position of the biomagnetic field component does not correspond to the position of the current source. Further, with only the isointegral map obtained from the electrocardiogram, because of the individual difference such as the position and size of the internal organs, making it difficult to accurately determine abnormality such as a disease by simply gathering from the isointegral map.
An object of the present invention is to provide biomagnetic field measuring method and apparatus which can grasp the whole state of a living body portion by using maps which are greatly reduced in number as compared to the maps required in the prior arts.
Another object of the present invention is to provide biomagnetic field measuring method and apparatus which can permit analysis of a magnetic field source by measuring a vertical component Bz of a biomagnetic field without increasing the number of detection coils.
According to the present invention, (1) a biomagnetic field measuring method comprises: a first step of measuring a temporal change of a component of a biomagnetic field generated from a living body by using a plurality of fluxmeters disposed externally of the living body and each including a superconducting quantum interference device (SQUID), the magnetic field component being in a first direction which is vertical to the surface of the living body; a second step of determining a temporal change of a value proportional to a root of square sum of change rates of the first-direction magnetic field component in second and third directions which cross the first direction; a third step of integrating the temporal change of the value obtained in the second step over a predetermined interval to determine an integral value, and a fourth step of displaying the integral value obtained in the third step.
According to the present invention, (2) a biomagnetic field measuring method comprises: a first step of measuring temporal changes of components of a biomagnetic field generated from a living body by using a plurality of fluxmeters disposed externally of the living body and each including a superconducting quantum interference device (SQUID), the magnetic field components being in first and second directions which are parallel to the surface of the living body; a second step of determining a temporal change of a value proportional to a root of square sum of the first-direction and second-direction magnetic field components; a third step of integrating the temporal change of the value obtained in the second step over a predetermined interval to determine an integral value; and a fourth step of displaying the integral value obtained in the third step.
Specifically, in the biomagnetic field measuring methods (1) and (2) as above, the above integral values are used through interpolation and extrapolation to display an isointegral map for connecting points at which the integral values in the above fourth step are equal to each other, the above third step of integrating the temporal change of the value obtained in the second step over a predetermined interval to determine the integral value is carried out over a plurality of predetermined intervals to determine a plurality of integral values, and computation for determining any of the ratio, the sum or the difference between the plurality of integral values is carried out. In the Cartesian coordinate system (x, y, z), the direction vertical to the body surface is defined as z axis, the first direction is defined as z direction, the second direction is defined as x direction and the third direction is defined as y direction. In the polar coordinate system (r, xcex8, xcfx86), the direction vertical to the body surface is defined as r axis, the first direction is defined as r direction, the second direction is defined as xcex8 direction and the third direction is defined as xcfx86 direction.
According to the present invention, (1) a biomagnetic field measuring apparatus for measuring biomagnetic field distribution comprises: a plurality of fluxmeters disposed externally of a living body and each including a superconducting quantum interference device (SQUID) for detecting, as a signal, a biomagnetic field generated from the living body; operation processing unit for performing the operation processing of the signal; and display unit for displaying a result of the operation processing. In the biomagnetic field measuring apparatus, the fluxmeters detect a temporal change of a component of a biomagnetic field, the magnetic field component being in a first direction which is vertical to the surface of the living body, the operation processing unit performs computation for determining a temporal change of a value proportional to a root of square sum of change rates of the first-direction magnetic component in second and third directions which cross the first direction and computation for integrating the temporal change of the value over a predetermined interval to determine an integral value, and the display unit displays the integral value.
According to the present invention. (2) in the above biomagnetic field measuring apparatus, the fluxmeters detect temporal changes of components of a biomagnetic field, the magnetic field components being in first and second directions which are parallel to the surface of the living body, the operation processing unit performs computation for determining a temporal change of a value proportional to a root of square sum of the first-direction and second-direction magnetic components and computation for integrating the temporal change of the value over a predetermined interval to determine an integral value, and the display unit displays the integral value.
Specifically, in the biomagnetic field measuring apparatus in (1) and (2) as above, an isointegral map for connecting points at which the integral values are equal to each other is obtained through interpolation and extrapolation and displayed on the display unit, and the operation processing unit carries out the computation of integrating the temporal change of the value over a predetermined interval to determine the integral value over a plurality of predetermined intervals to determine a plurality of integral values and computation for determining any of the ratio, the sum or the difference between the plurality of integral values, and the plurality of fluxmeters are disposed at equal intervals on the surface of the living body.
In the biomagnetic field measuring apparatus of the present invention, components of a magnetic field generated from the heart, that is a normal component a-d a tangential component which are respectively vertical and parallel to the chest surface can be displayed simultaneously. In the Cartesian coordinate system (x, y, z), when the direction vertical to the living body surface is assumed to be z axis, the first direction is defined as z direction, the second direction is defined as x direction and the third direction is defined as y direction. In the polar coordinate system (r, xcex8, xcfx86), when the direction vertical to the living body surface is assumed to be r axis, the first direction is defined as r direction, the second direction is defined as xcex8 direction and the third direction is defined as xcfx86 direction.
Essentially, in the present invention, when the direction vertical to the living body surface is assumed to be z axis of the Cartesian coordinate system (x, y, z) and the plane parallel to the living body surface is assumed to be (x, y) plane, a normal component Bz(x, y) of biomagnetic field vertical to the body surface is detected, and tangential components Bx and By of biomagnetic field parallel to the body surface are presumed from change rates of the normal component Bz in the x and y directions, respectively.
According to the present invention, without resort to detection coils for measuring the tangential components Bx and By, an isomagnetic map indicative of projection of current distribution upon the two-dimensional (x, y) plane can be obtained, a current source in the living body can be decided from a peak pattern in the isomagnetic map, and (x, y) coordinate positions of a plurality of current dipoles can be known.
The contents of the operation processing carried out by the operation processing unit (a computer such as a personal computer for collecting signals detected by a plurality of fluxmeters and applying the following operation processing to the collected signals or an electronic circuit in the form of hardware dedicated to the operation processing) will be described.
When a plurality of fluxmeters each including a superconducting quantum interference device (SQUID) are used to detect tangential components (parallel to the surface of a living body) Bx(x, y, t) and By(x, y, t) of a magnetic field generated from the living body at a position (x, y) on the body surface (where in the Cartesian coordinate system (x, y, z), the plane parallel to the body surface is assumed to be xy plane and the axis vertical to the body surface is assumed to be z), two-dimensional vector magnitude |Bxy(x, y)| (hereinafter, | | represents absolute value) is determined from a root of square sum of the tangential components Bx(x, y, t) and By(x, y, t) pursuant to equation (3).
|Bxy(x, y, t)|={square root over ( )}{(Bx(x, y, t))2 +(By(x, y, t))2}xe2x80x83xe2x80x83(3)
Subsequently, an integral value I1(x, y) of waveform |Bxy(x, y, t)| at each point (x, y) is obtained over a desired interval pursuant to equation (4), an isointegral map for connecting points at which the integral values I1(x, y) at respective points (x, y) are equal to each other is obtained through interpolation and extrapolation, and the isointegral map is displayed on the display screen.
I1(x, y)=∫|Bxy(x, y, t)|dtxe2x80x83xe2x80x83(4)
Hereinafter, presumption of the tangential components Bx and By from the measured magnetic field component Bz(x, y, t) vertical to the body surface will be described.
By taking advantage of the fact that the tangential component of biomagnetic field parallel to the body surface best reflects a current flowing through a portion immediately below the body surface and considering the relation between the current flow direction and the magnetic field direction, current distribution in the living body projected upon a two-dimensional plane parallel to the body surface can be surveyed by rotating a tangential vector (Bx, By) of the measured magnetic field counterclockwise through 90xc2x0. More particularly, where ex and ey represent unit vectors in x-axis and y-axis directions, a current vector J indicated by equation (5) can be determined from tangential components Bx and By at respective measuring points and can be expressed in terms of distribution (arrow map) of current vector fields at the respective measuring points (x, y).
J=xe2x88x92Byex+Bxeyxe2x80x83xe2x80x83(5)
On the other hand, when the normal component Bz of magnetic field vertical to the body surface, an arrow map using a current vector expressed by equation (1) is defined (the first prior art: H. Hosaka and D. Cohen (1976)).
J=(∂Bz/∂y)exxe2x88x92(∂Bz/∂x)eyxe2x80x83xe2x80x83(1)
Comparing equation (1) with equation (5), the present inventors have found the possibility that equations (6) and (7) stand, that is, the possibility that the tangential components Bx and By can be induced from the normal component Bz of the measured magnetic field and have studied in various ways. Results of studies will be described hereunder in greater detail.
Bx=xe2x88x92(∂Bz/∂x)xe2x80x83xe2x80x83(6)
By=xe2x88x92(∂Bz/∂y)xe2x80x83xe2x80x83(7)
FIG. 1 is a diagram useful for explaining modeling the generation of a magnetic field due to action of the heart (cardiac magnetic field) by a magnetic field generated from a current dipole in a horizontally layered conductor and analyzing the model. In FIG. 1, P designates a horizontally layered conductor having its surface on xy plane of the Cartesian coordinate system (x, y, z), Q designates moment of a current dipole existing at a position indicated by a position vector r0 (x0, y0, z0), and r(x, y, z) designates a position vector of a measuring point at which magnetic flux density B(r) (magnetic field) is measured. In the model shown in FIG. 1, a magnetic field B(r) generated outside the horizontally layered conductor P is formulated by Sarvas (literature: Phys. Med. Biol., Vol. 32, No. 1, pp.11-22 (1987)) and is expressed by equation (8).
B(r)={xcexc0/(4xcfx80K2)}{Qxc3x97axc2x7ez∇Kxe2x88x92Kezxc3x97Q}xe2x80x83xe2x80x83(8)
In equation (8), xcexc0 designates magnetic permeability of vacuum, ez designates a unit vector in z-axis direction, xe2x80x9cxxe2x80x9d designates vector product, xe2x80x9cxc2x7xe2x80x9d designates scalar product, and ∇ designates grad=(∂/∂x, ∂/∂y, ∂/∂z). Then, a is indicated by equation (9), a is indicated by equation (10), K is indicated by equation (11) and ∇K is indicated by equation (12). | | indicates absolute value.
a=r(x, y, z)xe2x88x92r0(x0, y0, z0)xe2x80x83xe2x80x83(9)
a=|a|xe2x80x83xe2x80x83(10)
K=a(a+axc2x7ez)xe2x80x83xe2x80x83(11)
∇K=(2+axe2x88x921axc2x7ez)a+aezxe2x80x83xe2x80x83(12)
Tangential components Bx and By of the B (r) given by equation (8) which are parallel to the horizontally layered conductor P and normal component Bz vertical to the horizontally layered conductor P are given by equations (13), (14) and (15), respectively.
Bx={xcexc0/(4xcfx80K2)}xc3x97[{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}(∇K)x+KQy]xe2x80x83xe2x80x83(13)
By={xcexc0/(4xcfx80K2)}xc3x97[{Qy(yxe2x88x92y0)xe2x88x92Qx(xxe2x88x92x0)}(∇K)x+KQx]xe2x80x83xe2x80x83(14)
Bz={xcexc0/(4xcfx80K2)}xc3x97[{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}(∇K)z]xe2x80x83xe2x80x83(15)
On the other hand, a differential in x direction of the normal component Bz indicated by equation (13) is expressed by equation (16).
∂Bz/∂x={xcexc0/(4xcfx80K2)}xc3x97[{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}
{xe2x88x922(∇K)z(∇K)x/Kxe2x88x92axe2x88x923(xxe2x88x92x0) (zxe2x88x92z0)2+axe2x88x921(xxe2x88x92x0)}xe2x88x92(∇K)zQy]xe2x80x83xe2x80x83(16)
Similarly, a differential in y direction of the normal component Bz is expressed by equation (17).
∂Bz/∂y=xe2x88x92{xcexc0/(4xcfx80K2)}xc3x97[{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}
{2(∇K)z(∇K)y/K+axe2x88x923(yxe2x88x92y0) (zxe2x88x92z0)2xe2x88x92axe2x88x921(yxe2x88x92y0)}+(∇K)zQx]xe2x80x83xe2x80x83(17)
In equations (16) and (17),
xcex1=(∇K)z/Kxe2x80x83xe2x80x83(18)
xcex2x=xe2x88x92axe2x88x923(xxe2x88x92x0) (zxe2x88x92z0)2+axe2x88x921(xxe2x88x92x0)xe2x80x83xe2x80x83(19)
xcex2y=xe2x88x92axe2x88x923(yxe2x88x92y0) (zxe2x88x92z0)2+axe2x88x921(yxe2x88x92y0)xe2x80x83xe2x80x83(20)
are put, equations (16) and (17) are reduced to equations (21) and (22).
∂Bz/∂x=xe2x88x92{xcexc0/(4xcfx80K2)}xc3x97[{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}{2xcex1(∇K)xxe2x88x92xcex2x}+xcex1KQy]xe2x80x83xe2x80x83(21)
∂Bz/∂y=xe2x88x92{xcexc0/(4xcfx80K2)}xc3x97[{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}{2xcex1(∇K)yxe2x88x92xcex2y}+xcex1KQx]xe2x80x83xe2x80x83(22)
For simplification, equations (13), (21), (14) and (22) are normalized by a common factor {xcexc0/(4xcfx80K2)} so as to be reduced to equations (13xe2x80x2), (21xe2x80x2), (14xe2x80x2) and (22xe2x80x2).
Bx=(∇K)x{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}+KQyxe2x80x83xe2x80x83(13xe2x80x2)
∂Bz/∂x=xe2x88x922xcex1(∇K)x{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}xe2x88x92xcex1KQy
+xcex2x{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)=xe2x88x922xcex1Bx+xcex1KQy+xcex2x{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}xe2x80x83xe2x80x83(21xe2x80x2)
xe2x80x83By=(∇K)y{Qy(yxe2x88x92y0)xe2x88x92Qx(xxe2x88x92x0)}+KQxxe2x80x83xe2x80x83(14xe2x80x2)
∂Bz/∂y=xe2x88x922xcex1(∇K)y{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}xe2x88x92xcex1KQx4
+xcex2y{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}=xe2x88x922xcex1By+xcex1KQx+xcex2y{Qx(yxe2x88x92y0)xe2x88x92Qy(xxe2x88x92x0)}xe2x80x83xe2x80x83(22xe2x80x2)
As will be seen from equations (13xe2x80x2) and (21xe2x80x2), the value of ∂Bz/∂x equals a value obtained by adding two additional terms to a term equal to xe2x88x922xcex1 times the tangential component Bx and as will be seen from equations (14xe2x80x2) and (22xe2x80x2), the value of ∂Bz/∂y equals a value obtained by adding two additional terms to a term equal to xe2x88x922xcex1 time the tangential component By.
When moment Q=(Qx, Qy, 0), where Qx=Qy=50[nAm], exists at a point r0(0, 0, xe2x88x92z0), where z0=0.05[m], inside the horizontally layered conductor P as shown in schematic positional relation of FIG. 2, Bx (equation (13)) is compared with xe2x88x92∂Bz/∂x (equation (16)). By substituting x0=y0=y=0 and Q0=0 to equations (13) and (16), equations (23) and (24) are obtained.
Bx(x, 0) =(xcexc0/(4xcfx80K2)}{xe2x88x92(∇K)xQyx+KQy}xe2x80x83xe2x80x83(23)
∂Bz(x, 0)/∂x={xcexc0/(4xcfx80K2)}{2xcex1(∇K)xQyxxe2x88x92xcex1KQyxe2x88x92xcex2xQyX}xe2x80x83xe2x80x83(24)
FIG. 3 shows Bx (equation (23)) and xe2x88x92∂Bz/∂x (equation (24)) on the horizontally layered conductor P in terms of relative magnetic field magnitude curves C1 and C2 which are normalized by maximum values of Bx and xe2x88x92∂Bz/∂x. More specifically, the curve C1 represents Bx(x, 0)/max|Bx(x, 0)| and the curve C2 represents {xe2x88x92∂Bz(x, 0)/∂x}/max|∂Bz(x, 0)/∂x|. As will be seen from FIG. 3, the distribution of each of the Bx and xe2x88x92∂Bz/∂x has a peak at the original (x=0) which is immediately above the existence of the current dipole, indicating that the maximum signals of both the Bx and xe2x88x92∂Bz/∂x can be detected when the measuring point is immediately above the point where the current dipole exists. The curve C2 has a sharper peak than the curve C1, indicating that the magnetic field distribution due to xe2x88x92∂Bz/∂x (equation (16)) has higher spatial resolution than the magnetic field distribution due to Bx (equation (13)).
Magnetic field magnitude curves C3, C4 and C5 depicted in FIG. 4 represent the first, second and third terms of xe2x88x92∂Bz(x, 0)/∂x, respectively. Gathering from the results shown in FIG. 4, the third term is negligible in relation to the first and second terms, so that the shape of xe2x88x92∂Bz(x, 0)/∂x can be deemed to be determined by the first and second terms and equation (24) can be approximated by equation (24xe2x80x2).
∂Bz(x, 0)/∂x=(xcexc0/(4xcfx80K2)}{2xcex1(∇K)xQyxxe2x88x92xcex1KQy)xe2x80x83xe2x80x83(24xe2x80x2)
FIG. 5 shows curves indicative of magnitude of relative magnetic field obtained by comparing the first term with the second term of each of the equations (13) and (16) after normalization. In FIG. 5, curve C6 represents {first term of Bx(x, 0)}/max|Bx(x, 0)|, that is, {xe2x88x92(∇K)xQyx}/max|Bx(x, 0)|, curve C7 represents {first term of xe2x88x92∂Bz(x, 0)/∂x}/max|∂Bz(x, 0)/∂x|, that is, {xe2x88x922xcex1(∇K)xQyx}/max|∂Bz(x, 0)/∂x|, curve C8 represents {second term of Bx(x, 0)}/max|Bx(x, 0)|, that is, {KQy}/max|Bx(x, 0)|, and curve C9 represents {second term of ∂Bz(x, 0)/∂x}/max|∂Bz(x, 0)/∂x|, that is, {xcex1KQy}/max|∂Bz(x, 0)/∂x|.
The results of FIG. 5 show that the distribution of each of the first and second terms of xe2x88x92∂Bz(x, 0)/∂x is sharper than the distribution of each of the first and second terms of Bx(x, 0) and the sharpness of the distribution is prescribed by xcex1=(∇K)z/K defined by equation (18).
In FIG. 6, magnetic field curve C10 represents xcex1=(∇K)z/K, magnetic field curve C11 represents xe2x88x92{first term of equation (24)}/{first term of equation (23)}, that is, 2xcex1(∇K)xQyx/(∇K)xQyx=2xcex1, and magnetic field curve C12 represents xe2x88x92{second term of equation (24)}/{second term of equation (23)}, that is, xcex1KQy/KQy=xcex1. As shown in FIG. 6, xcex1=(∇K)z/K (curve C10) has a peak point at the original where the current dipole exists, and the peak value is 2/(zxe2x88x92z0). The magnitude of xe2x88x92∂Bz(x, 0)/∂x differs from that of Bx(x, 0) by 2/(zxe2x88x92z0) at the peak point. The current dipole exists at a depth indicated by (zxe2x88x92z0). It is difficult to determine (zxe2x88x92z0) from practical measurement or magnetic field. By comparing equations (23) and (24xe2x80x2), equation (25) is obtained.                                           -                          ∂                                                B                  z                                ⁡                                  (                                      x                    ,                    0                                    )                                                              /                      ∂            x                          =                                                            {                                                      μ                    0                                    /                                      (                                          4                      ⁢                                              xe2x80x83                                            ⁢                      π                      ⁢                                              xe2x80x83                                            ⁢                                              K                        2                                                              )                                                  }                            ⁢                              xe2x80x83                            ⁢                              {                                                                            -                      2                                        ⁢                                          xe2x80x83                                        ⁢                                                                  α                        ⁡                                                  (                                                      ∇                                                          xe2x80x83                                                        ⁢                            K                                                    )                                                                    x                                        ⁢                                          Q                      y                                        ⁢                    x                                    +                                      α                    ⁢                                          xe2x80x83                                        ⁢                    K                    ⁢                                          xe2x80x83                                        ⁢                                          Q                      y                                                                                            )                    =                                    2              ⁢              α              ⁢                              xe2x80x83                            ⁢                                                B                  x                                ⁡                                  (                                      x                    ,                    0                                    )                                                      -                                          {                                                      μ                    0                                    /                                      (                                          4                      ⁢                                              xe2x80x83                                            ⁢                      π                      ⁢                                              xe2x80x83                                            ⁢                      K                                        )                                                  }                            ⁢                              xe2x80x83                            ⁢              α              ⁢                              xe2x80x83                            ⁢                              Q                y                                                                        (        25        )            
Namely, when the second term is smaller than the first term in equation (25), approximate equation (26) is deemed to stand.
xe2x88x92∂Bz(x, 0)/∂x=2xcex1Bx(x, 0)xe2x80x83xe2x80x83(26)
In generalization, when two additional terms other than xe2x88x922xcex1Bx is smaller than xe2x88x922xcex1Bx in equation (21xe2x80x2), approximate equation (27) is deemed to stand.
∂Bz/∂x=xe2x88x922xcex1Bxxe2x80x83xe2x80x83(27)
In the foregoing, the results of studies on the relation between xe2x88x92∂Bz/∂x and Bx are described but similarly, this holds true for the relation between xe2x88x92∂Bz/∂y and By and approximate equation (28) of equation (22xe2x80x2) is deemed to stand.
xe2x80x83∂Bz/∂y=xe2x88x922xcex1Byxe2x80x83xe2x80x83(28)
Hereinafter, the procedure for determining an isomagnetic map by presuming tangential components Bx and By from the measured normal component Bz on the assumption that Bx is proportional to xe2x88x92∂Bz/∂x and By is proportional to xe2x88x92∂Bz/∂y pursuant to equations (27) and (28) will be described in greater detail.
When a magnetic field component Bz(x, y, t) vertical to the surface of a living body is detected, change rate ∂Bz(x, y, t)/∂x in x direction of the Bz(x, y, t) and change rate ∂Bz(x, y, t)/∂y in y direction of the Bz(x, y, t) are determined and root S(x, y, t) of square sum of the change rates is determined as indicated by equation (33).
S(x, y, t)={square root over ( )}[{∂Bz(x, y, t)/∂x}2+{∂Bz(x, y, t)/∂y}2]xe2x80x83xe2x80x83(33)
Subsequently, a waveform St(t, x, y) at each point (x, y) is integrated over a desired interval to determine an integral value I2(x, y) pursuant to equation (34), and then an isointegral map for connecting points at which integral values I2(x, y) at the respective points (x, y) are equal to each other is obtained through interpolation and extrapolation and the isointegral map is displayed on the display screen.
xe2x80x83I2(x, y)=∫|St(x, y, t)|dtxe2x80x83xe2x80x83(34)
For example, when the heart is an object to be measured, intervals during which respective waves Q, R and S are generated, an interval during which a QRS wave (QRS complex) for generation of Q to S waves is generated and an interval during which a T wave is generated are used for the integration range in equations (4) and (34). Further, a plurality of integration ranges are taken in equations (4) and (34) to determine a plurality of integral values, computation for determining the sum, the difference or the ratio between the integral values is carried out, an isointegral map for connecting points at which the computation results have the same value is determined through interpolation and extrapolation, and the isointegral map is displayed on the display screen. For example, an interval T1 during which the QRS wave is generated is set as a first integration range and an interval T2 during which the T wave is generated is set as a second integration range, integral values I1,T1(x, y) and I2,T1(x, y) are determined for the interval T1 pursuant to equation (4) and integral values I1,T2(x, y) and I2,T2(x, y) are determined for the interval T2 pursuant to equation (34), and sum Isum(x, y) inclusive of isoweight, difference Idif(x, y) or ratio r(x, y) is determined between the integral values I1,T1(x, y) and and I1,T2(x, y) or between the integral values I2,T1(x, y) and I2,T2(x, y) pursuant to equations (35) and (36), equations (37) and (38) or equations (39) and (40).
Isum(x, y)=w1xc3x97I1,T1(x, y)+w2xc3x97I1,T2(x, y)xe2x80x83xe2x80x83(35)
Isum(x, y)=w1xc3x97I2,T1(x, y)+w2xc3x97I2,T2(x, y)xe2x80x83xe2x80x83(36)
Idif(x, y)=w2xc3x97I1,T2(x, y)xe2x88x92w1xc3x97I1,T1(x, y)xe2x80x83xe2x80x83(37)
Idif(x, y)=w2xc3x97I2,T2(x, y)xe2x88x92w1xc3x97I2,T1(x, y)xe2x80x83xe2x80x83(38)
r(x, y)=I1,T1(x, y)/I1,T2(x, y)xe2x80x83xe2x80x83(39)
r(x, y)=I2,T1(x, y)/I2,T2(x, y)xe2x80x83xe2x80x83(40)
The results of operations pursuant to equations (35) and (36), equations (37) and (38) and equations (39) and (40) suppress irregularities in isointegral map due to the individual difference and abnormality of living body function due to diseases can be detected.
With the isointegral map obtained in the present invention, states of all living body portions can be grasped by using the number of maps which is far smaller than the number of maps required in the prior arcs without analyzing biophenomena by the use of many maps, required in the prior arts, which indicate states of living body portions at respective time points. Since the peak position in the isointegral map obtained by using the tangential component or the normal component of biomagnetic field coincides with a portion in living body through which a large amount of current flows, it can be decided from the isointegral map which portion in living body a large amount of current flows within a desired time zone through. The biomagnetic field distribution differs greatly individual by individual but according to the present invention, the integral value over a desired interval obtained from a waveform representing a temporal change of a component in each direction of the biomagnetic field is used and therefore, more quantitative biomagnetic field distribution can be displayed by using a smaller number of maps and disease and abnormality of each individual can be grasped objectively and quantitatively.
Further, in the present invention, an isomagnetic map equivalent to the conventional isomagnetic map based on Bxy can be obtained by measuring only the normal component Bz without measuring tangential components Bx and By through vector measurement. With the conventional isomagnetic map obtained directly from the normal component Bz, a plurality of current sources are difficult to discriminate but in the isomagnetic map of the present invention, the peak pattern appears immediately above the current source as in the case of the conventional isomagnetic map based on Bxy, thus bringing about advantages that a plurality of current sources in the living body can be observed directly and the inverse problem of analyzing the position and size of the plurality of current sources can be solved with ease.
To summarize the present invention, reference is made to FIG. 7. More particularly, a biomagnetic field measuring apparatus of the present invention for measuring biomagnetic field distribution inside a shield room 1 has a plurality of fluxmeters each including a superconducting quantum interference device (SQUID) and operative to detect a biomagnetic field generated from a living body 2 in the form of a signal, an operation processing unit 8 for performing the operation processing of the signal, and a display unit for displaying results of the operation processing. The fluxmeters detect a temporal change of a vertical magnetic field component representing a component of the biomagnetic field in a first direction which is vertical to the surface of the living body, and the operation processing means performs computation for determining a temporal change of a value proportional to a root of square sum of change rates of the vertical magnetic field component in second and third directions which cross the first direction and computation for determining an integral value of the temporal change over a predetermined interval, and the display means displays the integral value. Since the quantitative biomagnetic field distribution is displayed by using a small number of maps, disease and abnormality of each individual can be grasped objectively and quantitatively.
Further, in the present invention, an isomagnetic map equivalent to the conventional isomagnetic map based on Bxy can be obtained by measuring only the normal component Bz without measuring the tangential components Bx and By through vector measurement and, by setting the number and position of peaks in a pattern of the obtained isomagnetic map to the initial condition, the inverse problem of analyzing the position and size of the current source in the living body can be solved with ease.